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ABSTRACT 

We use very deep near-infrared (NIR) imaging data obtained in MOIRCS Deep Survey (MODS) to 
investigate the evolution of the galaxy stellar mass function back to z ~ 3. The MODS data reach 
J = 24.2, H = 23.1, K = 23.1 (5a, Vega magnitude) over 103 arcmin 2 (wide) and J = 25.1, H = 23.7, 
K = 24.1 over 28 arcmin 2 (deep) in the GOODS-North region. The wide and very deep NIR data 
allow us to measure the number density of galaxies down to low stellar mass (1O 9 -1O 1O M0) even at 
high redshift with high statistical accuracy. The normalization of the mass function decreases with 
redshift and the integrated stellar mass density becomes ~ 8-18% of the local value at z ~ 2 and ~ 
4-9% at z ~ 3, which are consistent with results of previous studies in general fields. Furthermore, 
we found that the low-mass slope becomes steeper with redshift from a ~ —1.3 at z ~ 1 to a ~ —1.6 
at z ~ 3, and that the evolution of the number density of low-mass (10 9 -10 10 Mq) galaxies is weaker 
than that of M* (~10 11 M Q ) galaxies. This indicates that the contribution of low-mass galaxies to 
the total stellar mass density has been significant at high redshift. The steepening of the low-mass 
slope with redshift is opposite trend expected from the stellar mass dependence of the specific star 
formation rate reported in previous studies. The present result suggests that the hierarchical merging 
process overwhelmed the effect of the stellar mass growth by star formation and was very important 
for the stellar mass assembly of these galaxies at 1 < z < 3. 

Subject headings: galaxies: evolution — galaxies: high-redshift — infrared:galaxies 



1. INTRODUCTION 

Understanding when and how galaxies built up their 
stellar mass is one of the most important issues in 
observational cosmology and galaxy formation. Since 
the stellar mass of a galaxy is generally dominated 
by long-lived low-mass stars, it is considered as the 
integral of the past star formation rate. Then, the 
evolution of the cosmic stellar mass density provides 
a global picture of the past history of star formation 
in the universe. Recent near-infrared (NIR) surveys 
allow us to measure the average stellar mass density 
of th e universe directly from low to high redshift 
(e.g., IWilkins et al.l [20081 and references therein). At 
z < 1, many previous studies found the relatively mild 
evolution of the average stellar mass density, although 
there has been a significant evolution of the relative 
contributions from different morphological types or 
different S EP /color types to the tot a l stellar mass den- 
sity (e.g., iBrinchmann k. Ellis! |2000L iDrorv et all I2004L 
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I Bundv et al.l I2005L IPannella et al l I2006L iBorch et al 
20061 iFranceschini et all 120061 I Vergani et alj 1200 
lllbert et al.l l2009f K In contrast, at 1 < z < 3, a rapid 
evolution of the stellar mass density ha s been observed 
by s tu dies based on deep sur vevs (IDickinson et al.l 



20031 iFontana et all 120031 IRudnick et al.r~l200 
Fontana et al.ll2004lGlazebrook et al.ll2004L IDrorv et al 
20051 iGwvn fc Hartwicl I2005L IRudnick et all I2006L 
Fontana et al.ll2006l lArnouts et al.l 120071 . iPozzetti et al l 
2007L IP erez-G onzalcz et al .1 12008ft . and a significant 



fraction of the stellar mass in the present universe 
seems to have been formed at that epoch. This is 
qualitatively consistent with results that the average 
star formation rate density in the universe peaked 
around z ~ 1-3 (e.g . iMadau et alJ ll996t I Hopkins! l200l 
Hopk ins fc Beaco mT 120061) . although several studies 
pointed out a quantitative discrepancy between the evo- 
lution of the cosmic star formation rate density and the 
stella r mass density (IRudnick et al.ll2003l. IRudnick et al 
2006, Hopkins & Bca coml I2006L IP crcz-Gonzalez ct al 



2008, Wil kins et al.ll2008L see also lArnouts et al.ll2007f h 



Stellar mass function (SMF) of galaxies and its evolu- 
tion over cosmic time provide additional information on 
how stellar mass assembly of galaxies proceeded in dif- 
ferent mass ranges, while the evolution of the integrated 
stellar mass density enables us to investigate when stars 
formed in the universe as mentioned above. In paticu- 
lar, the evolution of the number density of massive (e.g., 
M s tar > 10 11 M Q ) galaxies has been extensively studied in 
order to address when massive galaxies assembled, which 
is considered as an important test for galaxy formation 
models (e.g.. iGlazebrook et al.ll2004 ICaputi et al.l l2006l. 
IConselice et al.l 120071 iBerta et all 120071 ). While the co- 
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moving number density of these massive galaxies shows 
no significant evolution at z < 1 , it decreases significantly 
with redshift at 1 < z < 3, although the observed number 
density of very massive galaxies at z ~ 2-3 seems to be 
still higher than th e predictions of the hierarchical galaxy 
formation models (iGlazebrook et aT1l2004l ICaputi et al.1 
[200llConselice et alJl2007ft~ In order to understand how 
the stellar mass assembly of massive galaxies proceeded, 
studies for the evolution of lower mass galaxies are also 
important, since many low-mass objects at high redshift 
are expected to merge over cosmic time to form mas- 
sive galaxies in hierarchical structure formation scenar- 
ios. The hierarchical merging process generally destroys 
small galaxies and builds massive galaxies, and therefore 
this process is expected to change the shape of the SMF, 
for example, to increase the characteristic mass of the 
Schcchter function and/or to flatten the low-mass slope 
with time. In order to address the shape of the SMF 
reliably, it is essential to have a large sample of galaxies 
over a wide range of stellar mass. 

To investigate the evolution of galaxies over a wide 
range of mass is also important in the context of 
the 'downsizing' scenario. In the local universe, mas- 
sive galaxies tend to have redder color and older 
stellar population, while most low-mass galaxies are 
young and actively star- forming (jKauffmann et al.ll2003_ 
iBrinchmann et afll2004 . The star formation histories 
of galaxies dep end st rongly on thei r stellar mass (e.g., 
iHeavens et aTT 12004 [Jimenez et all 12005( 1 . Therefore 
the histories of the stellar mass assembly of high and 
low-mass galaxies could also be different. For exam- 
ple, higher star formation rate relative to stellar mass 
(i.e., specific star formation rate) of low-mass galaxies 
is expected to lead to higher growth rate of their stel- 
lar mass than that of massive galaxies, if the effect of 
galaxy mergers is ignored. In fact, several studies re- 
ported that the most massive galaxies with M sta r ~ 
10 -10 Mg show milder evolution in number den- 
sity at 1 < z < 3 than galaxies with M sta r ~ 10 11 - 
10 115 M Q , i.e., the mass-dependen t evolution of the num- 
ber density of massive galaxies (IConselice et al.l I2007L 
iBertone et aT1l2007l iPerez-Gonzalez et al]|2008D . At z < 
1—1.3, it was also observed that lower mass galaxies with 
Mstar ~ 10 10 -10 1:L M^ evolve faster than massive galaxies 
(|Fontana et al.ll2006l iPozzetti et al.ll2007tl . 

However, it is difficult to extend the investigation of 
galaxies over a wide range of stellar mass up to z ~ 1-3, 
which is an important era in the histories of the stellar 
mass assembly of galaxies. This is because very deep 
NIR data are necessary to unbiasedly sample galaxies at 
high redshift down to low st ellar mass, e.g., 10 9 -10 10 M Q 
(|Kaiisawa fe Yamada 1 12006). NIR data sample the rest- 
frame optical to NIR for galaxies at z ~ 1-3. Since NIR 
luminosity reflects total stellar m ass of a galaxy relatively 
well (e.g., iBell fe de Jondl2001l ) and is not strongly af- 
fected by dust extinction, NIR data are more suitable 
for construction of stellar mass limited samples than 
shorter wavelength data. NIR data are also important in 
that they straddle the redshifted Balmer/4000A break of 
these high-z galaxies, which is correlated strongly with 
the stellar M/L ratio and enables us to estimate stellar 
mass with high accuracy. Since such deep NIR data have 
been very limited so far, only a few studies investigated 
the evolution of the shape of the SMF up to z ~ 1-3 



using ultra deep NIR data with relatively small sury ey 
areas (|Gwvn fc Hartwicldl2005l iMarchesini et alJl2008T >. 

In this paper, we use very deep and wi de NIR data ob- 
tained in MOIRCS Deep S urvey fMODS; lKaiisawa et all 
120061 [ichikawa et al. 2007^1 with publicly available multi- 
wavelength data i n the Great Observator ies Origins Deep 
Survey fGOODS. iGiavalisco et al] 12004 in order to in- 
vestigate the evolution of the SMF of galaxies back to 
z ~ 3. The wide and deep NIR data allow us to construct 
a statistically large galaxy sample over a wide range of 
stellar mass even at high redshift. Section 2 describes the 
data set and the procedures of source detection and pho- 
tometry We present methods for constructing the SMF 
and investigate its evolution in section 3. In section 4, 
we compare the results with previous studies and discuss 
their implications for galaxy formation and evolution. A 
summary is provided in section 5. 

We use a cosmology with H o =70 km s _1 Mpc -1 , fi m = 
0.3 and Oa = 0.7. The Vega-referred magnitude system 
is used throughout this paper. 

2. OBSERVATIONAL DATA AND SOURCE DETECTION 

In the MOIRCS Deep Survey, we obtained the very 
deep JHK s-bands imaging data in the field of GOODS- 
North with Multi-O bject InfraRed Camera and Spec- 
trogr aph (MOIRCS, llchikawa etall I2006L ISuzuki et aD 
2008) mounted on the Subaru Telescope. The obser- 
vations were carried out in the period from 2006 April 
to 2008 M ay. We also used a rchiv al MOIRCS data ob - 
tained by IWang et~aT1 (|2009H and iBundv et all (2009). 
Four MOIRCS pointings cover ~ 70% of the GOODS- 
North region (~ 103.3 arcmin 2 , hereafter referred as 
'wide' field) and the total exposure time is 6.3-9.1 hours 
in J-band, 2.5-4.3 hours in ii-band and 8.3-10.7 hours 
in Xs-band. One of the four pointings , which includes 
the H ubble Deep Field North (HDF-N, IWilliams et~aT~l 
11996ft . is also the ultra-deep field of the MODS (hereafter 
'deep' field), where the exposure time is 28.2 hours in J- 
band, 5.7 hours in i?-band and 28.0 hours in Ks-hand. 

The data were reduced using a purpose-made IRAF- 
based software package called MC SRED 8 . The reduction 
procedures were as described in iKaiisawa et al.l (|2006l ) 
except that a defringing proc ess was additionally ap- 
plied. A standard star P177-D (|Leggett et al.ll2006f l was 
used for flux calibration. In a photometric night, object 
frames for flux calibration were observed immediately af- 
ter the standard star at similar airmass in all 4 pointings 
in JH K s-bands. Before combining the object frames, 
all frames were scaled to the count level of these calibra- 
tion frames using relatively bright unsaturated sources. 
We discarded the frames where the count level was very 
low (less than 70% of the calibration frame) or the im- 
age quality was bad (FWHM of the PSF larger than 
1.2 arcsec for the wide field and 0.8 arcsec for the deep 
field). The combined data reach J = 24.2, H = 23.1 and 
K = 23.1 (5a, Vega magnitude) for the wide field and 
J = 25.1, H = 23.7 and K = 24.1 for the deep field. 
Further details of the observations, reduction and qual- 
ity of the data will be presented in a forthcoming paper 
(Kajisawa et al. in preparation). 

The source detection was performed in the iCs-band 
image using the SExtractor image analysis package 

8 http:/ /www. naoj.org/stafT/ichi/MCSRED/mcsred_e. html 
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Fig. 1. — Photometric redshift vs. spectroscopic redshift for galaxies with spectroscopic redshifts from the literature. Four panels show 
the different population synthesis models used as SED templates in the photometric redshift estimate. 



(jBertin fc Arnoutsl Il996l) . We adopted MAG.AUTO 
from the SExtractor as the total K-b&nd magnitudes of 
the detected objects. In this study, we use the magnitude 
limited samples with K < 23 and K < 24 for the wide 
and deep fields, respectively. The detection complete- 
ness for point sources is more than 90% at the if-band 
limits in the both samples and the false detection rate 
is expected to be less than ^ 1% from the simulation of 
the source detection on the inverse Xs-band image with 
the same detection parameters (Kajisawa et al. in prepa- 
ration). We detected 6402 and 3203 sources above the 
i^-band limits in the wide and deep fields, respectively. 

In order to measure the optical-to-MIR spectral energy 
distributions (SEDs) of the detected objects, we used the 
publicly available multi- wavelength data in the GOODS 
field namely KPNO/MOSAIC ([/-band, ICapak et all 
[200l . HST/ACS (g,^,i,;z-band s version 2.0 data, Gi- 
avalisco et al., in preparation, iGiavalisco et al.l l2004h 
and Spitzcr/IRAC (3.6/im, 4.5/im, 5.8/zm, DR1 and 
DR2, Dickinson et al., in preparation), as well as the 
MOIRCS J and f/-bands images. These multi-band im- 
ages were aligned to the -ftfs-band image. The ACS and 
MOIRCS images were convolved to match the data with 
the poorest seeing for the both wide (to ~ 0.6 arcsec 



FWHM) and deep (to ~ 0.5 arcsec FWHM) samples. 
For the color measurements, we used a fixed aperture 
size with a diameter of 1.2 (1.0 for the deep sample) arc- 
sec (~ 2x seeing FWHM) for the ACS and MOIRCS 
data. For MOSAIC ([/-band) and IRAC (3.6^m, 4.5/xm, 
5.8/im) data, we first performed the aperture photom- 
etry of the wide (deep) sample with aperture sizes of 
2.28(2.39), 2.89(3.07), 3.05(3.24), 3.80(4.02) arcsec, re- 
spectively, and then applied the aperture correction by 
using the light profiles of the B- and ifs-bands images 
smoothed to match the resolution of (7-band and IRAC 
images, respectively (see Kajisawa et al. in preparation 
for details). 

In the following analysis, we use objects which are de- 
tected above 2a level in more than three bands (Ks-b&nd 
and other two bands) because it is difficult to estimate 
the photometric redshift and stellar mass of those de- 
tected only in one or two bands. In the wide and deep 
-ftT-band magnitude limited samples, only 21/6402 and 
42/3203 were excluded by this criterion, respectively. As 
a result, the wide and deep samples consist of 6381 and 
3161 objects, respectively, and the total number of ob- 
jects used are 7563, in which 1979 objects with K < 23 
in the deep field are included in the both samples. 




Fig. 2. — Uncertainty of the estimated stellar mass as a function of stellar mass in each redshift bin. Red circles represent the wide 
sample and blue circles show the deep sample. Vertical dashed lines show the limiting stellar mass described in the text for the wide (red) 
and deep (blue) samples, respectively. Open squares represent the median values at each stellar mass for all sample. For objects without 
spectroscopic redshift, the photometric redshift error is taken into account in the estimate of stellar mass uncertainty (see text). 



3. ANALYSIS 

3.1. Estimate of redshift and stellar mass 

In order to estimate the photometric redshift and stel- 
lar mass of the sample galaxies, we performed SED 
fitting of the multi-band photometry described above 
{UBVizJHK, 3.6/xm, 4.5/^m, 5.8/im) with population 
synthesis models. We adopted the standard mini- 
mum x 2 method for the fitting procedure. The re- 
sulting best- fit redshift (i.e., photometric redshift) is 
used for objects without spectroscopic identifications 
and the best-fit stellar mass-to-luminosity (M/L) ra- 
tio is used to calculate the stellar mass. We adopted 
spectroscopic redsh i fts (if availabl e ) from the literature 
(ICohen et all [2000l ICohenl [20011 [Dawson et ail [2001, 



Wirth et all 12004 ICowie et all l200l iTreu et all 120051 
Chapman et al.l I2005L iReddv et al.ll2006al iBarger et al l 



20081) , including those from our NIR spectro scopic obser- 
vation of ~ 20 star-forming BzK galaxies (|Daddi et al.l 
12001 at z ~ 2 with Subaru/MOIRCS (Yoshikawa et al. 
in preparation). For these objects, the SED fitting was 
performed fixing the redshift to each spectroscopic value. 

In this study, t hree models, i.e., GALAXEV 
Bruzual fc Charfotl I2003 D. PEGASE version 2 



Fioc fc Rocca-Volmerangd Il997) and iMarastonl (|2005t ) 



model, were used as SED templates. Comparisons 
among the results with different SED template sets 
allow us to check the systematic effects of these models 
on the estimate of redshift and stellar mass. We 
also used a public p hotometric redshift code, EAZY 
(|Brammer et al.l 120081) for an independent check of the 
photometric re dshift. In all SED models, Salpeter IMF 
(|Salpeterjll955f ) with lower and upper mass limits of 0.1 
and 100 M is adopted for easy comparison among the 
results with the models and those in other studies. The 
details of the model templates from each population 
synthesis model are the following. 

(i) GALAX EV — We used the simila r template sets 
described in lKaiisawa fc Yamadal (|2005l ). Exponentially 
decaying star formation histories (SFHs) with the 
decaying timescale r ranging between 0. 1 and 20 Gyr 
are a ssumed. Calzetti extinction law (|Calzetti et al.l 
I2000T ) in the range of E(B - V) = 0.0-1.0 is adopted. 
Metallicity is changed from 1/50 to one solar metallicity. 

(ii) PEGASE2 — We used the same template sets as 
those described in lGrazian et al.1 (|2006l ). Eight exponen- 
tially decaying SFHs with various gas infalling and star 
formation rate timescales are assumed (t„= k = 0.1-5 Gyr 
and t sfr = 0.1-20 Gyr, see Table 4 in iGrazian et al.l 
2006 for details). The model includes the metallicity 
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Fig. 3. — Top and middle: Rest-frame U — V color distributions of wide (top) and deep (middle) samples for each redshift bin. Dashed 
lines represent the if-band magnitude limits (K = 23 for the wide sample and K = 24 for the deep one) at the central redshift of each bin. 
Dashed-dot curves show the 90 percentile of U — V color at each stellar mass. Bottom: Comparison of the 90 percentile of U — V color 
between the wide (solid line) and deep (short dashed line) samples. Long-dashed line and dashed-dot line show the i-sT-band magnitude 
limits for the wide and deep samples. 



evolution and dust extinction in a self-consistent way. 
Following I Grazian et al.l (|2006l ). we also added passively 
evolving templates with constant star formation rate 
(SFR) and several truncation ages, and dusty star- 
forming templates with a constant SFR and the Calzetti 
extinction with the range of E(B — V) = 0.5-1.1. 

(iii) Maraston (2005) - - The same SFHs and the 
same range of the Calzetti extinction as in the case 
of GALAXEV are assumed (r = 0.1-20 Gyr and 
E(B -V) = 0.0-1.0). Metallicity is changed from 0.5 
to 2.0 solar metallicity. This model increases emphasis 
on contributions of thermally pulsating AGB stars to 
the SEDs. In all these (i)-(iii) models, the model age 
is changed from 50 Myr to the age of the universe at 
the observed redshifts. We also added to the model 
SEDs the Lyman series abso rption produc ed by the 
intergalactic medium following iMadaul (|1995[ ). 

(iv) EAZY — We used the default parameter setting 
and the default template set of eazy_vl.0. The default 
template set consists of 5 principal component templates 
constructed from a large number of PEGASE models 
that cover a distributio n of the SFH of gal a xies in the 
semi-analytic model by iDe Lucia fc BlaizoH (|2007l ) . In 
addition to these principal components, a dusty star- 
forming SED is added as a complementary template. 
EAZY fits the multi-band photometry of the observed 



gala xies with a linear co mbination of these templates 
(see lBrammer et al.l [20081 for details). Since EAZY fits 
with the templates based on the principal components, 
we cannot derive the precise stellar M/L ratio from the 
fitting procedure. In order to estimate the stellar mass, 
therefore, we performed SED fitting with the (i)-(iii) 
models, fixing the redshift to the output from EAZY. 
The results with the different population synthesis 
models at the same redshifts also enable us to check 
the systematic effect of these models on the estimate of 
stellar M/L ratio independently. 

In Figure [TJ we compare the photometric redshifts 
with a sample of 2102 spectroscopic redshifts from 
the literature described above. The photometric red- 
shift accuracy and the fraction of the catastrophic fail- 
ure with (5z/(l+z spoc ) >0.5 for each SED model are 
(<5z/(l+z spcc ) = -0.009 ± 0.077, 4.0% outliers) for 
GALAXEV, (-0.002 ± 0.103, 4.0%) for PEGASE2, 
(-0.007 ± 0.082, 3.9%) for the Maraston model and 
(+0.003±0.098, 4.7%) for EAZY, respectively. Although 
the photometric redshift accuracy is relatively good in all 
cases, there are some systematic differences among the 
models especially at 1.5 < z < 3.0 in Figure [1] We 
check the effects of these differences of the photometric 
redshifts on the SMF in section 13.41 

Table [1] lists the number of objects in each redshift bin 
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Fig. 4. — Stellar mass limit as a function of redshift for the wide and deep samples. Three panels show the results with the different SED 
models. Solid lines show the wide sample and dashed lines show the deep one. 



TABLE 1 
Sample size in each redshift bin 



redshift model templates Wide a Deep 5 



z 


=0.5- 


-1.0 


GALAXEV 


1945 


(859) 


902 (311) 








PEGASE2 


1860 


(859) 


836 (311) 








Maraston 


2118 


(859) 


1023 (311) 








EAZY 


2141 


(859) 


991 (311) 


z 


=1.0- 


-1.5 


GALAXEV 


1426 


(353) 


635 (105) 








PEGASE2 


1319 


(353) 


608 (105) 








Maraston 


1403 


(353) 


573 (105) 








EAZY 


1141 


(353) 


464 (105) 


z 


=1.5- 


-2.5 


GALAXEV 


1306 


(209) 


666 (75) 








PEGASE2 


1380 


(209) 


677 (75) 








Maraston 


1015 


(209) 


488 (75) 








EAZY 


1334 


(209) 


774 (75) 


z 


=2.5- 


-3.5 


GALAXEV 


487 


(95) 


366 (57) 








PEGASE2 


507 


(95) 


370 (57) 








Maraston 


427 


(95) 


302 (57) 








EAZY 


494 


(95) 


356 (57) 



for different SED models. These redshift bins, of which 
width is sufficiently larger than the typical photometric 
redshift errors, were defined so as to include a reasonable 
number of galaxies in the bin for calculating the SMF. 
The comoving volumes are 8.5 x 10 4 Mpc 3 (2.5 x 10 4 
Mpc 3 for the deep field) at 0.5 < z < 1.0, 1.4 x 10 5 
Mpc 3 (3.7 x 10 4 Mpc 3 ) at 1.0 < z < 1.5, 3.4 x 10 5 Mpc 3 
(9.2 x 10 4 Mpc 3 ) at 1.5 < z < 2.5 and 3.4 x 10 5 Mpc 3 
(9.3 x 10 4 Mpc 3 ) at 2.5 < z < 3.5, respectively. The 
parenthesis in Table [1] represents the number of objects 
with spectroscopic redshift. The fraction of spectroscop- 
ically identified sources is relatively high, thanks to the 

extensive spectroscopic surveys in this field. 

As discussed in previous studie s (e.g.. iPapovich et al] 
120011 iKajisawa fc Yamadal I2005L IShaplev et all l2005f ). 
the uncertainty of the stellar mass is smaller than other 
parameters such as stellar age, star formation timescale, 
metallicity and dust extinction. While the broad-band 
SEDs are degenerated with respect to these parameters 
of the stellar populations (stellar age, star formation 
timescale, metallicity and dust extinction), the stellar 
M/L ratio (and therefore stellar mass) is much less af- 
fected by the degeneracy because the effects of these pa- 
rameters on the stellar M/L ratio tend to be canceled 
out with each other. Figure O shows the uncertainty of 
the stellar mass estimated from a \ 2 map in SED fitting 
with the GALAXEV model for each object. For sources 
without spectroscopic identifications, we also varied red- 



shift as a free parameter in the calculation of the x 2 map 
to take into account the photometric redshift error. The 
uncertainty increases with decreasing stellar mass and in- 
creasing redshift. The stellar mass errors become ~ 0.3- 
0.4 dex at the limiting stellar mass (the vertical dashed 
lines in the figure) described in the next subsection. 

3.2. Stellar mass-limited sample 

The iVband magnitude limited sample does not have 
a sharp limit in stellar mass even at a fixed redshift, be- 
cause the stellar M/L ratio at the observed iVband varies 
with different stellar populations. We used the rest-frame 
U — V color distribution as a function of stellar mass in 
each redshift bin to estimate the limiting stellar mass 
above which most of galaxies are expected to be brighter 
than the magnitude limits and detected in the fCs-band 
image. Since rest -frame U — V is corr e lated with stellar 
M/L ratio (e.g., iRudnick et~aTll2003L iMarchesini et all 
one can predict the mass dependence of stellar 
M/L and then estimate the effect of magnitude limit on 
the stellar mass distribution. 

Figure [3] shows the rest-frame U — V distribution in 
each redshift bin for the wide sample (top panels) and the 
deep (middle panels) sample. A dashed line in each panel 
represents the iVband magnitude limit (K = 23 for the 
wide sample and K — 24 for the deep one). All objects 
with stellar mass larger than this line (on the right side 
of the line in the figure) at each U — V value are brighter 
than the magnitude limit. In order to calculate the line, 
we used the stellar M/L ratio and the rest-frame color of 
the GALAXEV model with various star formation his- 
tories, dust extinction, metallicity. The maximum mass 
was selected from the possible range of the stellar mass 
for the models with K = 23 (or 24 for the deep sample) 
and each U — V color in order to depict the dashed line 
in the figure. Dashed-dot lines in Figure [3] represent the 
90 percentile of U — V color at each stellar mass. We 
adopted the point where the lines of the magnitude limit 
and 90 percentile of U — V color cross with each other as 
the limiting stellar mass. Above this limiting mass, more 
than 90% of objects are expected to be brighter than the 
X-band magnitude limit. 

The bottom panels in Figure [3] show the comparison 
of the 90 percentiles of U — V color for the wide and 
deep samples. The 90 percentiles for the both samples 
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Fig. 5. — Evolution of the stellar mass function of galaxies in the MODS field. From top to bottom, the panels show the results with 
the different population synthesis models. Circles and squares show the SMF calculated with the 1/V ma x formalism for the wide and deep 
samples. Open symbols indicate data points located below the limiting stellar mass, where the incompleteness could be significant. Solid 
symbols show data points above the limiting stellar mass. The results of the deep sample are plotted by shaded symbols at stellar mass 
where the wide sample is also above the limiting mass. Error bars are based on the Poisson statistics. The solid lines show the results 
calculated with the STY m ethod for the all ( wide and deep) samples. The best-fit Schechter parameters are also shown in each panel. For 
reference, the local SMF of Cole et al. (2001) is shown as the dashed line. 
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Fig. 6. — Evolution of stellar mass function for different SED models. For each case, only data points above the limiting stellar mass 
are plotted. At stellar mass where the wide and deep samples are above the limiti ng mass, th e data points for the wide sample are plotted 
(the same ones as solid symbols in Figure E}. Dashed line shows the local SMF of ICole et all pOOlT ). 



agree well with each other even near the limiting mass of 
the wide sample. This suggests that the incompleteness 
near the magnitude limit of the wide sample does not 
strongly affect the U — V color distribution, although 
some red galaxies might be missed on the left side of the 
long-dashed line. 

Figure [4] shows the calculated limiting stellar mass 
as a function of redshift for the wide and deep sam- 
ples. We also estimated these mass limits for the 
cases with PEGASE2 and Maraston models as well 
as the case with GALAXEV. As seen in Figure [3] 
and o ther previous studies (e.g., iKaiisawa fc Yamadaf 
20051 iKaiisawa &: Yamada I 120061 lLabbe et all 120051 
Taylor et al. 1 120091) . less massive galaxies tend to have 
bluer rest-frame color even at high redshift. Such a mass- 
dependent color distribution can be seen well above the 
X-band magnitude limit up to at least z ~ 2.5 in Fig- 
ure [3l Since the bluer color of low-mass galaxies indi- 
cates a lower M/L ratio, we can detect galaxies down to 
the relatively lower mass limit with high completeness 
compared with, for example, the mass limit based on 
the M/L ratio of the passively evolvin g models, which 
is used in other previous studies (e.g.. iDickinson et al.1 
120031 iFontana et al.l200l . On the other hand, our 'wide' 
field data are relatively shallow for galaxies at z > 3 (top 
right panel in Figured]) and the completeness is relatively 



low even at high mass, where galaxies tend to have red 
rest-frame colors (high M/L ratios). We use objects with 
the stellar mass larger than these mass limits at each red- 
shift to estimate and discuss the SMF in the following. 



3.3. Deriving the stellar mass function 

The SMF of galaxies was derived with the non- 
paramet ric 1/V ma x formalism and the parametric STY 
method (jSandage et al.l fl979h . Both methods are com- 
monly used to estimate the luminosity function and stel- 
lar mass function of galaxies. 

In the V max method, V max was calculated with the 
best-fit model SED template for each galaxy. For each 
best-fit SED, we estimated the isT-band apparent mag- 
nitude as a function of redshift, taking into account of 
both the luminosity distance and X-correction. Then we 
determined the maximum redshift, z max above which the 
object becomes fainter than the K-baiid magnitude limit 
(K = 23 for the wide sample or K = 24 for the deep 
sample). V max is a comoving volume integrated from 
the lower limit of each redshift bin to z max or the upper 
limit of the bin (the smaller of these two) . Then 1 /V max 
estimates were used to calculate the number density of 
galaxies in each mass bin. 

In the STY method, assuming the Schechter function 
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Fig. 7. — Best-fit Schechter parameters as a function of redshift. Different symbols represent the results with different population 
synthesis models . Data points for different SED models are plotted with small horizontal offsets for clarity. Those of the local SMF of 
ICole et al.l ( 2001) arc also shown for reference. 

with the stellar mass larger than Mi; m (z) in each red- 
shift bin. Both the wide and deep samples were used 
simultaneously in the maximum likelihood technique. 



form (|Schechterill976l ) for the SMF, we estimated best- 
fit values of the Schechter parameters (a, M* , <f>*). The 
limiting stellar mass Mi; m (z) (described in the previous 
subsection) for the redshift of each object was used to 
calculate the probability that the object has the observed 
stellar mass M as 



P = 



<j>(M) 



(1) 



Here 4>(M) is the SMF represented by the Schechter func- 
tion. We searched for the values of the Schechter param- 
eters (a, M*, </>*) which maximize the likelihood L= J|p, 
the products of the probability densities for the objects 



3.4. Evolution of the stellar mass function 

Figure [5] shows the stellar mass function (SMF) of 
galaxies in the different redshift bins for the different 
SED models. The results of the V max method and 
the best-fit Schechter function estimated with the STY 
method are plotted in each panel. Error bars are based 
on the Poisson statistics. Dashed lines show the lo- 
cal SMF derived from the 2dF and 2MASS surveys 
(|Cole et al.ll2001l) with the small correction for the "max- 
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TABLE 2 

Best-fit Schechter parameters obtained with the STY method 

SEP model rcdshift bin a log 10 M* (M Q ) log 10 <j>* (Mpc~ J ) 

GALAXEV z = 0.5-1.0 -1.26 fHf 11.33 f^S -2.79 

,= 1.0-1.5 -1.48 +°g 11.48 ±S:i| -3-40 

,= 1.5-2.5 -1.52 I";" 11-38 ±g;« -3.59 tg; 1 * 

Z = 2.5-3.5 1. 75 ■!;;,] 11-42^° lili^oj 

PEGASE2 , = 0.5-1.0 -1.21 fjfro 11.31 +008 -2.73 

2=1.0-1.5 -1.32 ±g : g* 11.36 ±g : ig -3.16±°-'i° 

,= 1.5-2.5 -1.45 ±°;g| 11.32 l*" -3-51 ±°i| 

z = 2.5-3.5 -1.59 tp.H 11-39 ±pjg -3.98 tgjg 

Maraston 2 = 0.5-1.0 -1.33 +003 11.43 toto -3.04 ™™ 

2 = 1.0-1.5 -1.42 ±°J 11.31 ±g:il -3.36 1*12 

2 = 1.5-2.5 -1.35 tool 1°- 98 to'io " 3 - 40 -oil 

2 = 2.5-3.5 -1.58 +g 11-03 ±g;^ -3.93 +gg 

EAZY + GALAXEV 2 = 0.5-1.0 -1.30 toot H-37 tooi -2.84 too's 

2=1.0-1.5 -1.35 ±°;g| 11-43 -3-28i° ii 

2=1.5-2.5 -1.58 ±g;g| 11.50 ±g;il -3-73 ±g;Jg 

2 = 2.5-3.5 -1.63 ±g;^| 11.33 ±g;|g -3.99 ±g;|l 

EAZY + PEGASE2 2 = 0.5-1.0 -1.30 +^ 11.45 +g ^ -2.87 

2=1.0-1.5 -1.29 1°^ 11.41 ±gj» -3.19 t^? 

2=1.5-2.5 -1.48 11.48 tg : » -3.62 1°;^ 

+0.15 11 Q7 +0.32 , nn +0.28 

ll.OI _ .22 -1.UU _ 39 

EAZY + Maraston 2 = 0.5-1.0 -1.34 to 03 H-45 toot -3.06 

2=1.0-1.5 -1.38 ±£°j 11.39 ±g|| - 3 - 42 j:i4 

z — 1.5 2.5 -1.52 q 11.12 Q -^-^ _ 3.55 q -yj 

2 = 2.5-3.5 -1.59 ta 1 ^ 10.98 t°o' 2 2o " 3 - 84 io.38 



imum age" method as described in lFontana et al.l <|2004f h 
In Figure O we plot the combined wide and deep com- 
plete data (same as the solid symbols in Figure [5]) for the 
different SED models with different symbols in the same 
panel. 

Figures [5] and [6] show that the SMFs obtained from dif- 
ferent samples and SED models are in good agreement, 
although there are some systematic differences among 
the SMFs. The number densities for the deep sample 
are systematically larger (by ~ 0.1 dex) than those for 
the wide sample at 0.5 < z < 1.0 for all SED models. 
Our deep field is centered at the HDF-N field where the 
extensive spectroscopic surveys reve aled the large scale 
structures at z = 0-85 and z = 1.02 (jCohen et al.ll2000L 
iWirth et al.l 120041 ). These filaments or clumps around 
the HDF-N could cause the slightly larger number den- 
sity in our deep field, and such differences can be con- 
sidered as the possible field-to-field variance. Slighter 
excess in the deep field can be also seen at 2.5 < z < 3.5 
in most SED-model cases, and this may also be due to 
large scale structures. On the other hand, the system- 
atic differences among the different SED models in Fig- 
ure [S] seem to be larger. The number density of galaxies 
for the Maraston model is systematically smaller by ~ 
0.15-0.2 dex (up to ~ 0.5 dex at > 10 n M Q ) than those 
for the GALAXEV and PEGASE2 models especially at 
z > 1.5. Since the differences are also seen between the 
EAZY+Maraston and EAZY+GALAXEV/PEGASE2, 
where the same redshifts of EAZY are used, the differ- 
ences of the estimated stellar M/L ratio among the SED 
models cause the sys tematic differences of the number 
density. For example, Mar aston et al.l (|2006h performed 
the broad-band SED fitting of relatively young (~ 0.2- 



2 Gyr stellar age) galaxies at 1.4 < z < 2.7 with the 
Maraston and GALAXEV models and reported that the 
Maraston model gives systematically younger age and 
lower stellar mass (~ 60%) than the GALAXEV model. 
Such difference of the estimated stellar mass seems to 
explain the differences seen in Figure [6] If the fraction 
of young galaxies becomes larger at high redshift, the 
larger differences in the SMFs at z > 1.5 could be also 
explained because the contribution of TP-AGB stars is 
expected t o be significant in the relatively young ages (~ 
0.2-2 Gvr. !Marastor]l2005L iMaraston et al.ll2006h . 

We can see general evolutionary features in Figures 
[5] and [6] in spite of the field variance and the system- 
atic differences among the SED models mentioned above. 
We note, first, that the overall number density decreases 
gradually with redshift in all cases. While the number 
density of galaxies at 0.5 < z < 1.0 is similar with that in 
the local universe, the number density at 2.5 < z < 3.5 is 
about an order of magnitude smaller than the local value. 
Figure [7] and Table [2] show the best- fit Schechter param- 
eters (a, M*, (j)*) estimated with the STY method. The 
rcdshift evolution of the overall number density can be 
seen as the decrease of the normalization of the SMF <f>* . 
4>* decreases down to ~ 50% of that in the local universe 
at z ~ 0.75, - 16% at z ~ 1.25, ~ 9% at z ~ 2 and - 3% 
at z ~ 3. Similar evolution of the S MF is also seen in pre- 
vious studies of gener al fiel ds fe.g.. iFontana et al.l l 2006l 
iPerez-Gonzalez et al.ll2008L iMarchesini et al.ll2008l ). 

Second, we found the mass-dependent evolution of the 
SMF. The evolution of the number density of low-mass 
galaxies with M sta r ~ 1O 9 -1O 1O M is smaller than that of 
massive galaxies with stellar mass ~ 10 11 M Q . While the 
number density of galaxies with ~ 10 1:l Mq at 2.5 < z < 
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Fig. 8. — Evolution of the Schechter parameters in M*-c* plane for the different SED models. Crosses show the best-fit values determined 
with the STY method, lcr (solid) and 2<r (dashed) error contours are also shown. 



3.5 is smaller by a factor of ~ 20 than the local value, 
that of galaxies with ~ 5 x 1O 9 M is smaller by only 
a factor of ~ 6 at the same redshift. The shape of the 
SMF at 0.5 < z < 1.0 is similar with that in the local 
universe, and it becomes steeper with redshift at z > 1. 
This can be seen in Figure [7] and Table [2] as a steepening 
of the low-mass slope a with redshift. a decreases with 
redshift gradually from a = -1.29±0.03(±0.04) at 0.5 < 
z < 1.0 to a = -1.48±0.06(±0.07) at 1.5 < z < 2.5 and 
a = -1.62 ± 0.14(±0.06) at 2.5 < z < 3.5. The quoted 
errors are statistical errors estimated from the maximum 
likelihood method. The values in parenthesis show the 
uncertainty due to the different SED models. 

The uncertainty of the Schechter parameters becomes 
larger with redshift, especially at 2.5 < z < 3.5, because 
of the larger limiting mass at higher redshift and of the 
small number of galaxies even at relatively high mass due 
to the evolution of the overall number density mentioned 
above. Nonetheless, since the uncertainty at z < 2.5 is 
rather small by virtue of our deep and wide NIR data, 
the evolution of the low-mass slope of the SMF between 
;$ z < 3 is found to be significant. Figure [5] shows 
the best-fit Schechter parameters and their uncertainty 
for the different SED models in the M*-a plane, which 
represents the evolution of the shape of the SMF. In all 
cases, the evolution of a is significant, although there is 
degeneracy between M* and a especially at 2.5 < z < 
3.5, where we can reach only to the relatively high stellar 
mass. 

On the other hand, the characteristic mass M* shows 
no significant evolution except for the results with the 
Maraston model. The M* values at 0.5 < z < 3.5 are 
similar with or slightly larger than those in the local 
universe. No significant evolution for M* is also seen 



in previous studies dFontana et al.l 120061 iPozzetti et al.l 
120071 iMarchesini et al.ll2008fl . For the Maraston model, 
M* becomes smaller by a factor of ~ 2-2.5 at z > 1.5. 
As in the above discussion of the overall number den- 
sity, this can be explained by the systematically lower 
stellar M/L ratio of the Maraston model because the 
same result is also seen when the same photometric red- 
shifts (EAZY) are used. While passively evolving galax- 
ies dominate the massive end (> 10 1:L Mr?)) of the SMF 
at z < . 1 (e.g.. I Juneau et alj [20051 iBorch et alj 120061 
IVergarri et ail 12008. Il bert et alj l2009f). many massive 
star- forming (i.e., rel atively young) galaxies have been 
found at z ~ 2 (e.g.. iDaddi et aTfl2007l iPapovich et al.1 



[2CM iBorvs et all 120051 . iShaplev et alJl2004h . It is pos- 
sible that TP-AGB stars contribute to the SED of mas- 
sive galaxies significantly only at z > 1.5. Therefore, the 
Maraston model would give systematically lower stellar 
mass for these massive galaxies. 

3.5. Possible biases for the evolution of the low-mass 

slope 

Here we investigate possible biases for the results in 
the previous subsection, in particular, systematic effects 
which could cause the steepening of the low-mass slope 
at high redshift. 

The larger limiting mass at higher redshift causes the 
degeneracy between M* and a as mentioned in the pre- 
vious section. Since low-mass galaxies near the limiting 
mass tend to be faint in each band, the photometric er- 
rors are larger, which results in the large uncertainty in 
their photometric redshift. The large errors of the pho- 
tometric redshifts of faint objects might lead to the sys- 
tematic increase of the low-mass galaxies at high redshift 
because the redshift distribution of the X-selected sam- 
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Fig. 9. — Monte Carlo simulation for the effects of photometric redshift uncertainty on the shape (M* and a) of stellar mass function 
with the GALAXEV model. Large square in each redshift bin shows the observed values at 0.5 < z < 1.0 and these are assumed not to 
evolve with redshift in the simulation. Small circles show the results of 200 simulations (see text for details). Observed M*— a values (cross) 
and la (solid) and 2<r (dashed) error contours are also shown for each redshift bin. 



pie has a peak around z ~ 1 and a tail to higher redshift. 
In order to evaluate the effect on the low-mass slope, we 
performed a Monte Carlo simulation, assuming no evo- 
lution of the shape of the SMF (i.e., constant M* and a). 
At first, we constructed mock catalogs with the observed 
M* and a at 0.5 < z < 1.0. For the normalization of 
the SMF, we assumed 0* evolving as (f>*(z) oc (1 + z)~ 2 
so that the redshift distribution of the mock sample 
is consistent with our observation. Even if we assume 
4>*{z) oc (1 + z)^ 1 or <f)*(z) oc (1 + z)~ 3 , the results shown 
in the following do not change significantly. Stellar mass 
and redshift of mock objects were randomly selected from 
the ranges of 10 8 M Q < M sta r < 10 12 M and < z < 6, 
using the probability distribution estimated from the as- 
sumed SMF and the corresponding comoving volume of 
our survey at each redshift. For each mock object, we 
randomly extract an object from the observed sample 
with similar mass and redshift (allowing duplicate), and 
adopted its observed multi-band photometry. The ob- 
served multi-band photometry was extracted from the 
catalog that contains all sources detected on the Ks-band 
image (including K > 23 or K > 24 objects) in order to 



take into account the scattering of objects fainter than 
the magnitude limits. Then we added random offsets to 
the multi-band photometry according to the measured 
photometric errors, and adopted the mock object if the 
resulting if-band magnitude of the object was brighter 
than the magnitude limits (K < 23 for the wide sam- 
ple and K < 24 for the deep sample). We repeated this 
procedure and made the mock catalogs with the same 
sample sizes as the observed wide and deep samples. 
The same SED fitting procedure as for the observed one 
was performed in order to estimate the photometric red- 
shift and stellar mass of the mock objects. For objects 
with spectroscopic identification, redshifts are fixed to 
the spectroscopic values. We performed 200 simulations 
and calculated the best-fit Schechter parameters in each 
simulation with the STY method. 

Figure [9] shows the results of the simulation in the case 
with GALAXEV (the results for the other models are 
similar). A relatively large scatter of a is seen in the 
highest redshift bin, which probably reflects the large 
degeneracy between M* and a due to the large limiting 
mass. Furthermore, the simulated a distributes around 
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Fig. 10. — Effect of the catastrophic failure of the photometric redshift on the SMF in the Monte Carlo simulation. Upward and downward 
triangles show the fractional increase and decrease of the number of galaxies as a function of stellar mass in each redshift bin due to the 
objects whose photometric redshift was changed catastrophically (<5z/(l + z) > 0.5) by the random offsets of the multi-band photometry. 
Upward (downward) triangles represent the fraction of the objects which enter into (drop out from) the redshift bin due to the catastrophic 
failure. Data points and error bars represent the median value and 68 percentile interval of the 200 simulations. Vertical lines show the 
limiting stellar mass for the wide (solid) and deep (dashed) samples. 



a systematically steeper (by ~ 0.1-0.15) value than the 
assumed one at 2.5 < z < 3.5, while the simulated val- 
ues tend to be slightly flatter (by ~ 0.1) in lower redshift 
bins. However, since the observed evolution of a is much 
stronger than the systematic effects in Figure [9j the ob- 
served steepening of the low-mass slope at high redshift 
is significant, especially at 1 < z < 2.5. 

In the simulation, the random offsets added to the 
multi-band photometry and the recalculation of the pho- 
tometric redshift include the effect of the catastrophic 
failure. We extracted the mock objects whose photomet- 
ric redshift was changed catastrophically (8z/(l + z) > 
0.5) by the random offsets and checked the effect of these 
objects on the resulting SMF. Figure fTOl shows the frac- 
tional increase and decrease of the number of galaxies 
due to the catastrophic failure of the photometric red- 
shift as a function of stellar mass in each redshift bin. 
At 0.5 < z < 1.5, there is ~ 10-20% decrease near the 
limiting stellar mass, while there is only negligible effect 
of the catastrophic failure at M star > 10 10 Mq. Most 
contamination from z < 0.5 or z > 1.5 occurs only at 
the stellar mass lower than limiting mass and it is not 
plotted in the figure. About 10-20% decrease near the 



limiting stellar mass might cause a slightly flatter low- 
mass slope seen in Figure but the effect is relatively 
small (~0.1 dex decrease in the number density). At 
z > 1.5, the effect of the catastrophic failure is similar or 
even smaller than that at low redshift. Furthermore, the 
contamination from lower redshift and the dropout from 
the redshift bin tend to be canceled out with each other, 
which results in the negligible effect on the SMF. 

In Figure [11] we also show the results of the same sim- 
ulation on the a~4>* plane to investigate the effect on the 
normalization. The systematic effect on the normaliza- 
tion of the SMF due to the photometric redshift errors 
seems to be relatively small, although the degeneracy be- 
tween the parameters makes direct comparison difficult. 

Next, we discuss a possibility of the over-deblending of 
faint objects with relatively low S/N ratio near the detec- 
tion limit. Although the -ftf-band, where we performed 
the source detection, corresponds to the rest-frame B- 
band even at Z = 3.5, the morphological K-correction 
could enhance the over-deblending at high redshift be- 
cause galaxies tend to show the patchy appearance in 
shorter wavelengths due to the dominance of young 
stars and the dust extinction (e.g.. lKuchinski et al.ll200l[ 
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Fig. 11. — Effect of photometric redshift uncertainty on the 
normalization (<j>*) of stellar mass function. Large squares and 
small circles are the same as those in Figure [9] but different colors 
represent different redshift bins. The normalization is assumed to 
evolve with redshift as 4>*{ z ) oc (1 + z)~ 2 (see text). 



iRawat et a.1.1 l2009h . Figure [T2l shows the fraction of the 
objects with the deblending flag by SExtractor as a func- 
tion of stellar mass in each redshift bin for the wide and 
deep samples. We cannot see the mass nor redshift de- 
pendence of the fraction of the deblended objects. We 
conclude that the effect of deblending does not cause the 
evolution of the low-mass slope. 

Finally, we estimated the stellar mass of galaxies in 
the highest redshift bin with the GALAXEV templates 
with two-component (old and young) star formation his- 
tories. If the old stellar population is hidden by re- 
cent star formation, the stellar mass could be underes- 
timated especially for relatively blue low-mass galaxies 
at hi gh redshift, where the S/N ratio tends to be low 
(e.g. JPapovich et al.ll200ll iDrorv et al.ll2005l ). We used 
the exponentially decaying star formation models (young 
component) with an old population component. For the 
young population, free parameters are stellar age, star 
formation timescale r, color excess, metallicity (same as 
for one-component star- formation history). For the old 
component, we limited the star formation timescale and 
stellar age to shorter and older values than the young 
population, respectively, and assumed no dust extinc- 
tion. Figure Q2] shows the comparison of the stellar 
masses estimated with one- and two-component mod- 
els for galaxies at 2.5 < z < 3.5. No significant sys- 
tematic difference of the stellar mass can be seen. The 
scatter is consistent with the uncertainty of the stellar 
mass of these galaxies shown in Figure [2j although the 
stellar mass with the two-component model is slightly 
larger for a small fraction of galaxies at low-mass end. 
The long wavelength data with Spitzer/IRAC, which 
sample the rest- frame NIR region even at high redshift, 
could make th e syst ematic uncer tainty relatively small 
(|Fontana et al.ll2006t lElsner et al.ll2008[ ). 

4. DISCUSSION 
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Fig. 12. — Fraction of deblended objects as a function of stellar 
mass in each redshift bin for the wide (top) and deep (bottom) 
samples. Solid line shows all i\-selected galaxies in each redshift 
bin and shaded histogram represent deblended objects. Dashed- 
dot line shows the fraction of the deblended objects. Vertical long- 
dashed line shows the limiting stellar mass. The results with the 
GALAXEV model are shown. 

We have studied the evolution of the SMF of galax- 
ies back to z ~ 3 using the NIR Subaru/MOIRCS data 
obtained in the MOIRCS Deep Survey and the publicly 
available multi-wavelength GOODS data. Our deep and 
wide NIR data allowed us to estimate the number den- 
sity of galaxies down to low stellar mass (~ 1O 9 -1O 1O M0) 
even at z ~ 2-3 with high accuracy. Main results in the 
previous section are 1) the decrease of the overall number 
density of galaxies with redshift, and 2) the steepening 
of the low- mass slope of the SMF at high redshift. In 
this section, we first compare the results with previous 
studies and then discuss the implications for galaxy for- 
mation and evolution in the following. 

4.1. Comparison with previous studies 

In this subsection, we compare our results in the 
MODS field with those in other general fields. 

Figure [TJ] shows a comparison of the integrated stel- 
lar mass density as a function of redshift. We calcu- 
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Fig. 13. — Comparison between the stellar masses estimated 
with one-component SFH model and two-component (old and 
young) model for galaxies at 2.5 < z < 3.5. The result with the 
GALAXEV model is shown. 



lated the stellar mass density by integrating the best- 
fit Schechter function over 10 8 -10 13 M in each redshift 
bin. Different symbols represent the re sults with the 
differe nt SE P models. Com p ilation s of iWilkins et al.l 
(|2008f) and IMarchesini et al.l (|2008f) are also plotted. 
All data points are the results with the Salpeter IMF. 
The stellar mass density in the MODS field decreases 
with redshift as seen in other previous studies. When 
compared with those in the local universe, the stellar 
mass density is ~ 53-72% of the local value at 0.5 < 
z < 1.0, - 22-34% at 1.0 < z < 1.5, - 8-18% at 
1.5 < z < 2.5 and ~ 4-9% at 2.5 < z < 3.5. For the 
GALAXEV mod el, which is widely used in other previ- 
ous studies (e.g.jDrorv et"al . 2005, Fontana et al. 120061. 
Pozzetti et al.ll2007. Eis ner et al. 2008, Marchesin i et all 



2008), our results are consistent with those of previous 
studies at each redshift (Figure [14|) . A slightly large 



value at 0.5 < z < 1.0 is probably due to the large scale 
structures around the HDF-N mentioned above. As dis- 
cussed in Section 13.41 the Maraston model gives system- 
atically lower stellar mass density, which is about three- 
quarters of those with the other two models (GALAXEV 
and PEGASE2) at 0.5 < z < 1.0, about two-thirds at 
1.0 < z < 1.5 and about a half at z > 1.5. 

In Figure [T5] we compare the SMF of galaxies in the 
MODS field with those in other studies. For simplic- 
ity, we plotted the result with the GALAXEV model. 
At the massive end (> 10 11 M Q ), our result is consis- 
tent with those in other general fields, although the 
uncertainty of our data and the variance among the 
different surveys are relatively large because of the 
small number of observed massive galaxies and prob- 
ably their strong clustering fe.g.. | van Dokkum et all 
120061 iConselice et "all [20071 IStringer et all I2009D . The 
agreement of the SMF at the massive end is consis- 
tent with the comparison of the integrated stellar mass 
density because galaxies around M* (~ 10 11 M W ) dom- 
inate the stellar mass density (e.g., iGwvn fc Hartwickl 
120051 see also next subsection). No significant evo- 
lution of M* seen in the MODS field at 0.5 < z < 
3.5 is also consistent with results of previous studies 



On the other hand, at the lower mass of the SMF, there 
is some variance in different fields including the MODS 
field in spite of the relatively small statistical errors. In 
addition to the differences in the normalization, which 
are similar offsets with those at the massive end, there is 
differences in the low-mass slope especially at high red- 
shift. One possible reason for this is the different depths 
of the surveys and therefore the different limiting stellar 
masses. The SMF in the MODS field shows a slight up- 
turn around ^10 10 M Q and a steeper slope at < 10 10 M Q , 
although this is marginal in the highest redshift bin. This 
indication can be seen in the results with all SED mod- 
els (Figure O and does not depend on the systematic 
error of the photometric redshift and stellar M/L ra- 
tio. Similar upturns in the S MF are also seen in field 
galaxies in the local universe (jBaldrv et all 12001) and 
at intermediate and high redshift dFontana et alj 120061 
IPozzetti et "all 120071 lElsner et aLll2008fh Therefore the 
relatively shallow survey with the limiting stellar mass 
larger than ~ 10 10 M Q could miss the upturn and result 
in the systematically flatter low-mass slope even if the 
completeness to the limit is sufficiently high or the in- 
completeness correction is properly done. For example, 
the SMFs from the different surveys show consistent low- 
mass slopes in the 0.5 < z < 1.0 bin in Figure [TBI where 
all the surveys reach below 1O 1O M . At 1.0 < Z < 1.5, 
several surveys reach only to ^1O 1O M . While the slope 
at 1O 1O -1O 1:L M0 is similar with those in other surveys, 
the SMF in the MODS field shows the steeper slope at 
10 9 -10 10 Mq than that at larger stellar mass. This could 
be the case in the higher redshift bins. 

At z > 1.5, on the other hand, lElsner et all (|2008l ) 
and IMarchesini et all (|2008l ) reach below 10 10 M Q and 
show the fl atter low-mass slopes than those in the present 
study. In lElsner et al.l (|2008l ). the low-mass slope is 
fixed to the weighted mean value a ~ —1.36, which 
is constrained mainly by galax ies at relatively low red - 
shift (z < 1, their Figure 7). IMarchesini et all (|2008l ) 
discussed the effects of various random and system- 
atic uncertainties on the estimate of the Schechter pa- 
rameters of the SMF and suggested that the uncer- 
tainty of the low-mass slope in their analysis is rela- 
tively large especially at high redshift due to the statis- 
tical error and the systematic errors of the different SED 
models. The difference in the low-mass slope between 
IMarchesini et alj (j2008) and the present study could be 
attributed to field-to-field variance. In the local universe, 
several studies suggest that the sha pe of the SMF de- 
pends on the environment of galaxies (|Balogh et aT1l200lL 
IBaldrv et aIll2006L IBaldrv et al.ll2008fl. and this could b e 
th e case at higher r e dshift (e.g.. iScodeggio et~aT1 [20091. 
In IMarchesini et alj (|2008l ) , the low-mass slope at high 
redshift is constrained mainly by the FIRES data, whose 
survey area is about a fifth of our study. If the shape 
of the SMF depends on the environment even at high 
redshift, the low-mass slope estimated from the small 
area surveys could be affected by the field-to-field vari- 
ance. A larger area survey with a similar or fainter lim- 
iting magnitude in the NIR wavelength is needed to con- 
clude whether there exists the field-to-field variance of 
the shape of the SMF at z ~ 2-3 or not. 
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FlG. 14. — Evolution of the stellar mass density (integrated over 10 8 Mq <M s tar < 10 13 Mq) as a function of redshift for the different SED 
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4.2. Mass-dependent evolution of galaxies 

The decrease of the overall number density and the 
steepening of the low-mass slope with redshift seen in 
the evolution of the SMF of galaxies in the MODS field 
suggest the possibility of the mass-dependent number 
density evolution. Figure [16] shows the evolution of 
the number density of galaxies in different stellar mass 
ranges for the different SED models. Relatively low- 
mass galaxies with the stellar mass of 10 9 ' 5 -10 10 ' 5 M Q 
tend to show a weaker evolution of the number den- 
sity especially at z > 1.5 than M* (~ 10 11 M Q ) galaxies, 
although the trend is relatively marginal in the results 
with the GALAXEV/PEGASE2 models. The relative 
number density of low-mass galaxies had been larger at 
high redshift than that in the local universe. In other 
words, the number density of galaxies around M* might 
have more rapidly increased at z > 1.5 than that of 
low-mass galaxies. On the other hand, several studies 
based on wide area but shallow surveys suggest that for 
M star > 10 n M Q , more massive (e.g., > 1O 115 -1O 12 M ) 
galaxies show weaker evolution in the number density a t 
z > 1.5 (|Berta et al.l 120071 iPerez-Gonzalez et al J 12008( 1 . 
Although the survey volume of the present study is not 
large enough to constrain strongly the number density 
of very massive galaxies with 10 115 -10 12 Mq, galaxies 
around 10 11 M Q might have been formed effectively at 
1 < z < 3 compared with higher and lower mass ga laxies 
(e.g., iFontana et al.ll2004l iFranceschini et al"]l2006( l. 

Figure [17] shows the contribution of galaxies in dif- 



ferent stellar mass ranges to the cosmic stellar mass 
density. While the stellar mass density is dominated 
by massive galaxies with M sta r > 1O 1O 5 M0 in the lo- 
cal universe, the contribution of lower mass galaxies 
increases with redshif t and is significant at z ~ 3. 
iReddv fc Steidell (|2009( 1 reported a relatively steep faint- 
end slope (a ~ ~l-7) of the luminosity function for Ly- 
man Break Galaxies (LBGs) at 1.9 < z < 3.4. They 
estimated the stellar mass of faint LBGs with the log- 
lin ear relation between SFR and stellar mass derived 
in ISawicki et al.1 (|2007), and suggested that those faint 
L BGs, which typically h ave M star < 10 10 M Q (Figure 12 
in lReddy fc Ste idcl 2009) , show a significant contribution 
to the total stellar mass density. Since most low-mass 
galaxies at z ~ 2-3 have rather b lue rest-frame colors 
(typi c ally rest U — V <0, Figure p iKajisawa fc Yamadal 
120051 IKajisawa fc Yamada l l2006Fand satisfy the LBG 
criteria ijLabbe et al.ll2005lv ! the steep low-mass slope of 
the SMF at high redsh ift in our study seems t o be consis- 
tent with the result in lReddv fc Steidell (|2009( 1. Note that 
the integrated stellar mass density in the MODS field is 
consistent with those of other previous studies as men- 
tioned in t he previous subsec t ion, b ut about a half of the 
estimate in lReddv fc Steidell (]2009l ). Then, although the 
contribution of low-mass galaxies to the total stellar mass 
density in the MODS field is larger than those estimated 
in previous studies, the stellar mass density is still smaller 
than that ex pected from the evolution of the cosmic SFR 
density (e.g.. lWilkinslt r a^[2008irReddv fc Steidel[2009ll . 
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Fig. 15. — Comparison between stellar mass functions in this study and other surveys. The GALAXEV model is used for the di- 
rect comparison with previous studies. Circles, squares and solid line are the same as those in Figure l5l For other deep surveys, the 
best-fit Schechter func tions from the literature HFontana et al.H2009 , lElsner et"aT]|200a fPozzetti et"alll2007l [Perez-Gonzalez et alJI2008L 
IMarchesini et al.|[2008T) are plotted over the observed mass range. 



4.3. The origin of the evolution of the low-mass slope 

What is the origin of the evolution of the low-mass 
slope ? In hierarchical structure formation scenarios, 
low-mass objects continue to grow by merging over cos- 
mic time and result in more massive objects at later 
epochs. Therefore one can expect the steeper low- mass 
slope at the high er redshift (e.g., iKhochfar et al.1 [20071 . 
iRvan et aI1l2007l) . Figure [H] compares the observed SMF 
in the MODS field with the predic tions of the publicly 
avai lable semi- ana l ytic m odels by iBower et all (|2006l ) 
and iBertone et all (|2007t ). These two models are based 
on the Millennium Simulation of the g rowth of dark 
matt er structure in a ACDM cosmology (Springel et al. 
2005) but are made by the different procedures in many 
aspects such as, e.g., the construction of halo merger 
trees, the implementation of the star formation, the 
feedback by supernova explosions and AGN activities. 
We convolved these model predictions with the Gaus- 
sian representing the uncertainty of the stellar mass es- 
timate shown in Figure [5J The median value in Fig- 
ure [2] is chosen as the width of the Gaussian depending 
on galaxy's stellar mass and redshift. In the models, 
the low-mass slope similarly steepens with redshift in 



Figure |TH] However, it should be noted that there are 
discrepancies between the observation and the models 
in the detailed shape of the SMF and that the evolu- 
tion of the overall number density of model galaxies is 
significantly weaker (i.e., larger number density at high 
redshift) than that of the present observa tion especially 
at low mass as in previous studi e s (e.g., Fontan a~et al.l 
[20061 IKitzbichler fc White! [20071 IBertone et afl [2001 . 
On the other hand, several studies found that semi- 
analytic models tend to predict the smaller number den- 
sity of galaxies at the massive e nd (M star > 10 11 ' 5 - 
10 1175 M f 7 1 ) than ob servations (e.g JConselice et al.|[2007l 
IBertone et al.|[2007t ). A similar indication is also seen in 
Figure [T51 although our data cannot strongly constrain 
the number density of these very massive galaxies. 

Since the stellar mass of a galaxy is generally domi- 
nated by long-lived low-mass stars, it can be generally 
considered as the integral of the past star formation 
rate and expected to increase monotonically with time. 
Therefore, for example, only suppression of star forma- 
tion in low-mass galaxies cannot cause the flattening of 
the low-mass slope of the SMF with time, while it could 
lead to the evolution of the faint-end slope of the rest- 
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Fig. 16. — Evolution of the number density of galaxies in different stellar mass ranges. Three panels represent the results with the 
different SED models. Circles and squares show the wide and deep samples, respectively. Different colors of the symbols represent different 
ranges of the stellar mass of galaxies. Only stellar mass ranges which are above the limiting stellar mass in Figure [4] at each redshift are 
shown. 
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Fig. 17. — Contributions of different stellar mass ranges to the total stellar mass density as a function of redshift. Symbols are the same 
as in Figure [TBI Diamonds show the stellar mass density integrated over 1O 9 M <M star < 1O 12 M at each redshift. 



frame UV luminosity function, which directly reflects 
the st ar formation rate, as discussed in iReddv fe Steidell 
(2009). On the other hand, the stellar mass growth by 
star formation activities in each galaxy could change the 
slope of the SMF. There are several observational studies 
that suggest that low-mass galaxies have higher star for- 
mation rate relative to their stellar mass (i.e., SFR/M st ar, 
specific star format ion rate) than mass i ve galaxies even 
at z ~ 2-3 (e.g. iFeulner et ~all 120051 IPapovich et all 
120061 IReddv et all l2006bl see "also ISantini et alj l2009f )7 
These results for the mass dependence of SFR/M sta r sug- 
gest that the expected growth rate of stellar mass by star 
formation activities in a galaxy is higher in lower mass 
galaxies. If we assume the observed trends of the low- 
mass slope and mass dependence of SFR/M sta r continue 
down to lower mass, a net increase of the number of 
galaxies in a mass bin due to star formation activities 
(inflow from lower mass bins — outflow into higher mass 
bins) is expected to be larger at lower stellar mass, which 
leads to the steepening of the low-mass slope with time, 
i.e., op posite to that seen in F igures [5] and [6l For ex- 
ample, iDrorv fe Alvare j (|2008f ) calculated the expected 
change of the SMF due to star formation only from the 
average star formation rate at each stellar mass of g alax- 
ies in the FORS Deep Field (jFeulner et all l2005h and 



showed that the number density of lower mass galaxies 
is expected to increase more rapidly with time especially 
a t z ^5 3. By subtracting the expected contribution due 
to star for mation from the observe d redshift evolution 
of the SMF. IDrorv fc Alvared (|2008f ) evaluated the effect 
of the hierarchical merging on the evolution of the SMF 
(i.e., destroying small galaxies and building more massive 
ones). In their analysis, at M sta r ^ 10 10 M Q , lower mass 
galaxies tend to be consumed in mergers at higher rate 
to counterbalance the effect of the mass-dependent star 
formation rate. Sin ce no evolution of the lo w-mass end 
slope is assumed in IDrorv fc Alvare3 (|2008h . the steep- 
ening with redshift indi cates the effect of mergi ng larger 
than that estimated in IDrorv fc Alvarei (|2008l ). There- 
fore the evolution of the SMF of galaxies in the MODS 
field suggests that the hierarchical merging process was 
important for the stellar mass assembly of these low-mass 
galaxies at z ~ 1-3, when the cosmic stellar mass density 
had been increasing rapidly. The relatively high fraction 
of peculiar/irregular morphol ogy in low-mass galaxies at 
z > 1 in the NIR wav elength (|Kajisawa fc Yamadall2005l 
IConselice et al.l f2005f ) might also imply intense merging 
activity of these galaxies. 

In this context, the fact that low-mass galaxies 
tend to have relatively younger stellar age and higher 
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Fi g. 18, — Compariso n of the observed stellar mass function wit h the predictions by t he semi-analytic models. Short-dashed line refers to 
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as in Figure [6] 



SFR/M s tar than massive galaxie s at low to intermediate 
redshift (i.e., 'down sizing', e.g.. iKauffmann et al.1 120031 . 
iHeavens et al.ll2004D might not necessarily mean relative 
scarcity of low-mass galaxies at high redshift. A sig- 
nificant fraction of old stars in relatively massive (e.g., 
~ 10 11 M Q ) galaxies at low redshift could be formed at 
high redshift in actively star-forming low-mass galaxies 
which had coalesced into more massive objects through 
the mergers. One of the means to address this is to in- 
vestigate the distribution of the SFR and mean stellar 
age of galaxies as a function of stellar mass down to low 
mass at high redshift. The combination of the distribu- 
tion of the SFR and mean stellar age as a function of 
stellar mass and the SMF shown in this study allows us 
to measure the contributions of galaxies in a certain mass 
range to the cosmic star formation rate density and stel- 
lar mass density and their age distribution as a function 
of redshift. We will be able to use these quantities in 
order to estimate what the mass range of high-z galaxies 
could have assembled into more massive galaxies at later 
epochs. 

5. SUMMARY 

In this paper, we have investigated the stellar mass 
function of galaxies at 0.5 < z < 3.5, using the very 



deep and wide NIR imaging data obtained in MOIRCS 
Deep Survey and the multi-wavelength public data from 
the GOODS. The MODS data reach K - 23 over ~ 103 
arcmin 2 and K ~ 24 over ~ 28 arcmin 2 . We constructed 
a large sample of galaxies down to ~ 1O 9 -1O 1O M0 up to 
z ~3. 

Our main results are as follows: 

• The normalization of the SMF decreases with red- 
shift gradually and the integrated stellar mass den- 
sity becomes ~ 8-18% of the local value at z ~ 2 
and ~ 4-9% at z ~ 3. The evolution of the stellar 
mass density estimated with the GALAXEV model 
is consistent with those in previous surveys. 

• The characteristic mass M* of the best-fit Schechter 
function shows no significant evolution. With the 
Maraston model, however, M* becomes smaller by 
a factor of ~ 2-2.5 at z > 1.5. This may be because 
many massive galaxies are in active star-forming 
phase at z > 1.5 and TP-AGB stars would con- 
tribute to the SED of these massive galaxies signif- 
icantly only at z > 1.5. 

• The low-mass slope of the SMF becomes steeper 
with redshift gradually from a = —1.29 ± 
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0.03(±0.04) at 0.5 < z < 1.0 to a = -1.48 ± 
0.06(±0.07) at 1.5 < z < 2.5 and a = -1.62 ± 
0.14(±0.06) at 2.5 < z < 3.5. The evolution of the 
number density of low-mass (10 9 -10 10 M Q ) galaxies 
is weaker than that of M* (~ lO n M ) galaxies. 
The contribution of these low-mass galaxies to the 
cosmic stellar mass density increases with redshift, 
and becomes significant at z ~ 3. 

• We noted a marginal upturn around ~ 10 10 M Q and 
a steeper slope at < 10 10 M Q in the SMF. Shallow 
data, which does not reach to < 10 10 M Q , could 
lead to underestimation of the low-mass slope. 

• The steepening of the low-mass slope with redshift 
could be explained as a result of the hierarchical 
merging process. Since the mass dependence of the 
SFR/M s tar distribution seen in previous studies is 
expected to lead to the steepening of the low-mass 
slope with time, the opposite trend found in this 
study suggests that the hierarchical merging pro- 
cess was very important for the stellar mass assem- 
bly of relatively low-mass galaxies at 1 < z < 3. 
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